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We develop an exact solution to the problem of one dimensional chiral bosons interacting via an 

s-wavc Feshbach resonance. This problem is integrable, being the quantum analog of a classical two- 
wave model solved by the inverse scattering method thirty years ago. Its solution describes one or 
two branches of dressed chiral right moving molecules depending on the chemical potential (particle 
density). We also briefly discuss the possibility of experimental realization of such a system. 



INTRODUCTION 



Hamiltonian takes the following form 



The advances of the last decade in the techniques of 
atomic physics allowed to realize a variety of exactly solv- 
able many-body models experimentally. In particular, 
the Tonks- Girardeau [1] gas, the dilute gas of repulsive 
one-dimensional bosons has been realized [2]. More gen- 
erally, the Lieb-Liniger model [3] of interacting bosons in 
one dimensions can now be studied experimentally. 

In the field of cold atomic gases it is customary to 
use Feshbach resonances to control interactions between 
the atoms. On some level, Feshbach resonances can be 
thought of as simply a tool to change the interaction 
strength. But on a deeper level, they are a way to convert 
pairs of atoms into molecules whose binding energy can 
be controlled. When confined to one dimensions, such a 
system is then described by the Hamiltonian [4] 
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Here o^, o are the creation and annihilation operators 

of atoms, w , b are those of molecules, rria, rrih arc their 
respective masses, and eo and g are two parameters con- 
trolling the resonance. Throughout the paper we take 
both atoms and molecules to be bosons. 

The model described by Eq. (1) is unlikely to be inte- 
grable even classically [5] and cannot be solved exactly. A 
particularly straightforward argument against quantum 
integrability involves calculating the amplitude of three 
atom scattering with incoming momenta pi, P2, P3 into 
outgoing momenta p'l , P2 , p-^ , distinct from any permuta- 
tion of pi, p2, P3, in the first nonvanishing Born approx- 
imation. This amplitude can be verified to be nonzero, 
while integrability would require it to be zero [6]. 

The bosonization techniques of developed in the con- 
text of Luttinger liquid theory can be applied to under- 
stand Eq. (1) [7]. An approximate technique based on the 
ideas of the asymptotic Bethe ansatz can also be used [8] . 

Yet there exist an exactly solvable model closely re- 
lated to Eq. (1). It is the model of chiral atoms and 
molecules interacting via a Feshbach resonance. Its 
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Eq. (2) describes the atoms which move in one direction 
with the velocity u, independent of their wave vector, 
and molecules which also move in one direction with the 
velocity v, also independent of their wave vector. The 
case of u = t; is degenerate, as we will see below. In what 
follows, it will be assumed that ti 7^ v. 

The classical version of this model is called the 2-wave 
model in the literature. It is known to be integrable and 
describes various phenomena, for example, in non-linear 
optics [9]. 

In this paper we demonstrate that this problem is also 
integrable quantum mechanically. We do that by em- 
ploying the coordinate Bethe Ansatz and developing the 
exact solution of the problem defined by Eq. (2). 

A closely related integrable model, called the 3-wave 
model [9] (which differs from Eq. (2) by having 3 chi- 
ral fields two of which can fuse into the third one) was 
studied and demonstrated to be integrable quantum me- 
chanically in Ref. [10]. This model is also interesting for 
the applications in quantum optics and atomic systems. 
Ref. [10] did not work out the finite density behavior for 
the 3-wave model, which is something we do here for the 
2-wave model. Studying this is thus an interesting direc- 
tion of further work. 

One can argue that it is not entirely straightforward 
to realize the Hamiltonian Eq. (2) using real atoms and 
molecules. However, this Hamiltonian can be thought 
of as an approximation to a true problem of atoms and 
molecules propagating in one dimension, given by Eq. (1), 
if we can restrict our attention to atoms and molecules 
whose momentum is close to a specially chosen momen- 
tum pp (for atoms) and 2pp (for molecules). Indeed, in 
that case their spectrum is linear, as in 
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We can now interpret pp/mas the velocity v the atoms in 
Eq. (2). One can immediately see a particular difficulty 
with this interpretation: the mass of the molecules has 

to be twice that of an atom, and so are their "Fermi" 
momenta. Thus the velocity of an atom and a molecule 
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must be equal. To justify the assumption that u ^ v, 
we may have to place the atoms and the molecules on 
an optical lattice where the molecular and atomic matrix 
elements are distinct from each other. Then the effective 
masses of atoms and molecules no longer have to be the 
same. 

The interaction in the Hamiltonian Eq. (2), controlled 
by the cubic Feshbach term, occurs at a point in space. 
In reality, however, the interactions between the atoms 
which lead to this term have a finite range which we 
denote Tq. Remembering that Vq is not zero is important 
in what follows. 



II. THE COORDINATE BETHE ANSATZ 

We begin by constructing a few body eigenstatcs of 
the Hamiltonian Eq. (2) and then proceed to generalize 
them to many body states made possible by the inte- 
grability of this problem. Then we impose the periodic 
boundary conditions to find the excitation spectrum of 
the system. In all steps we closely follow the standard 
techniques originally developed for the solution of the 
Lieb-Liniger model [3] as explained in Ref. [11]. 



A. Single atom state 

First of all, we observe that a single atom represents an 
exact eigenstate of the Hamiltonian Eq. (2). This state 
can be written as 



dx e'''''a^{x) |0) 
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where |0) is the vacuum, and its energy is given hy E = 
uk. Indeed, all terms in the Hamiltonian besides the very 
first one annihilate this state. 



B. Two atom states 

Now we seek the two-atom states in the following gen- 
eral form 
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Acting on it by the Hamiltonian Eq. (2) wc find the fol- 
lowing first-quantized Schrodinger equation 
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This equation has two classes of solutions. 

On the one hand, two atoms separated by some dis- 
tance will forever move at equal velocity u regardless of 



their momenta. Thus they will never interact. This is 
reflected in the existence of an exact eingenstate of the 
Hamiltonian given by ij){xi,X2) = f {\xi — X2\), (j) = Q, 
or 
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Here f{x) is an arbitrary function of its argument such 
that /(O) = 0. The energy of this state is simply E = 
u{ki + k2)- We call this a two-atom state. 

On the other hand, two atoms located at the same 
point will forever bo located at the same point since 
they move at an equal velocity. They hybridize with a 
molecule to form a dressed molecule state. It is given by 
ip{xi,x2) = AS{xi - X2)e^^^'^, (f){y) = Me^^^, equivalent 
to the following form 



/ 



dx e 



ikx 



Aa^^{x)+Mb\x)\ |0), 



(8) 



where A and M are the atom-mlecule amplitudes. Sub- 
stituting these into Eq. (6) gives 



ukA + -^M = EA, 
vkM + V2g6{0)A = EM. 



(9) 



Notice the appearance of 6{0). Treated naively, this could 
be interpreted as infinity. However, for the regularized 
model where the interactions happen at some length scale 
To, this term should be interpreted as (5(0) = 1/ro, which 
is what we do in what follows. We now find from Eq. (9) 
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Thus the dressed molecules come in two different species, 
labelled by the superscripts -|- and — . These two solution 
are easiest to parametrize if one introduces new variables 
X{k) which we call repidities (for reasons which will be- 
come clear later). We define those as in 
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Assuming, without the loss of generality, that v > u, X 
is positive if the plus sign is chosen in front of the square 
root and it is negative if the minus sign is chosen, thus 
in terms of A, the two species of dressed molecules are 
straightforward to distinguish (and so, the superscript ± 
is not necessary for A). In terms of these, one can write 
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as well as 



= eo(A), eo(A) = - (14) 

Here positive A correspond to the choice of + in the su- 
perscript and negative A to the choice of — . 

For completeness, we also note that k can be expressed 
in terms of A via 
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As seen from this relation, the energy of the dressed 
molecule is no longer a linear function of k, thus the 



dressed molecules have a nontrivial dispersion. As a re- 
sult, their velocity deo/dk depends on their wave vector 
k, and unlike atoms, the molecules can catch up with 
each other and undergo scattering, as we will see in the 
next subsection. 



C. Many-atomic states 

As happens in all integrable systems, more general 
states can be reduced to a combination of two body 
states. The most general eigenstate of the Hamiltonian 
Eq. (2) state can be written in the following form 
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where ^(a;) is the theta-function of its argument or 6{x) = 
1 if X > and 6{x) = if a; < 0. This state represents 
dressed molecules and M free atoms scattering off each 
other. An and M„ are the amplitudes for the dressed 
molecules corresponding to the rapidities A„, which in 
turn depend on the momenta fc„ (and on the species of 
the molecule, or on whether A„ is positive or negative). 
/ (Ij/a ~ 2//3I) is an arbitrary function of all possible differ- 
ences of the atomic coordinates, such that it vanishes if 
any two atomic coordinates coincide. Snm is an S-matrix 
representing the scattering of two dressed molecules, n 
and rn, off each other. Finally 5"™ is an S-matrix for 
scattering between an atom and a dressed molecule, with 
the only index n used to emphasize that it depends only 
on the rapidity of the molecule and not of the atom. In- 
deed, these S matrices take the following form 

_ i (A„ - A,„) + 1 _ 2iA„ - 1 

>->nm ~~ ■ / \ \ \ iJ'^n ii' \'- ' ) 

* \ ^ Kn) ^ 1 Z?,A„ + i 

This form of the scattering matrices justifies the term 
"rapidity" for the parameter A. It is now a matter of 

a straightforward algebra to check that Eq. (16) is an 
eigenstate of the Hamiltonian Eq. (2) with the energy 

N M 

i; = ^eo(A„)+M^p„. (18) 

n=\ n=l 

D. Bethe equations 

Construction of the exact cigenstates is but the first 
step towards exact solution of an integrable problem us- 



(16) 

I 

ing the coordinate Bethe ansatz. The next step is the 

imposition of the appropriate boundary conditions, the 
determination of the ground state energy and of the en- 
ergy of the excitations above the ground states. As usual, 

we impose the periodic boimdary conditions to arrive at 
the Bethe equations (here L is the system size) 

-i^n^i'=i'^"^"n^r=i- (19) 

gam (igpgnds on \j only and is pj independent. Thus it 
is always possible to choose pj in such a way that the sec- 
ond equation in (19) is satisfied. Then the first equation 

reduces to 

kjL + Y^ Oji - Me]"" = 2Tmj , (20) 

where iOji = In Sji, iO]"" = In 5?™. Following Ref. [11] 
it is straightforward to prove that the solution to these 
equations are unique and real, and all nj are distinct. 
This last claim is the consequence of the ID "Pauli" 
principle (at work here, as well as in the standard Licb- 
Linigcr model), which says that no two A can be the 
same, or the wave function Eq. (16) vanishes if \j = A/j 
for j ^k as, can be checked directly. 

As a next step, we take rij to be a continuous variable 
n{j), with Xj and kj becoming functions of n. This gives 

Lk{n) + 9{X{n) - X{1)) - M6''""(A(n)) = 27rn. (21) 

Here e{x) = \n[{ix + l)/{ix - 1)] /i and 6l'""(a;) = 
ln[{2ix-l)/{2ix + l)] /i 
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Finally, wc differentiate with respect to A(n), introduce 
the function 



P = 



1 dn 
LdX 



(22) 



playing the role of the density of A, and replace summa- 
tion by integration to arrive at 

P(A) 27ry^(A-M)2 + l~27rdA + 27rLl 



4A2' 
(23) 



E. Yang- Yang equation 

There are two ways to make further progress in the 
determination of the ground state energy and the exci- 
tation spectrum of the system. One follows explicit con- 
structions of the excitations, by exciting a state with a 
particular rapidity A, while shifting the rest of the rapidi- 
ties to accommodate the Bethe equations Eq. (20). The 
second is by studying the excitations at finite tempera- 
ture and then taking the limit T ^ 0. Both methods are 
described in Ref [11]. It is technically easier to use the 
second approach. Although this method is well known, 
we go over it briefly in the particular case of interest here. 

First we note that in a general state rij takes values 
in some subset of all possible integer numbers. We in- 
troduce pp as the density of A among the values of n 
which are taken ("occupied"), and ph as the density of A 
where these values are unoccupied, with pt = Pp + Ph (see 
Ref. [11] for the discussion on how this is done). Then 
we find 



I dk M 4 
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Next we construct the energy, the entropy, and the par- 
ticle number of such configuration, given by 

E = L J dXpp{X)eo{X), 
S = L j dX{ptlnpt- pplnpp- phlnph) , 
N = L J dXpiX). (25) 

Then we minimize the thermodynamic potential = 
E — TS — hN (T is the temperature, and h is chemical 
potential) with respect to pp, while remembering that the 
variation of Spp is related to 5pt by 
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Following standard methods [11], we introduce 
Ph 
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e(A) plays the role of the excitation spectrum of the sys- 
tem. It satisfies, as a result of the minimization of fl, 



<X) + lJdp 



(A - /i)2 + 
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Finally, we take the limit of zero temperature, T — >■ 0. 
This gives the following equation 



1 f ^. 2e(/x) 

Je(u.)<0 



e(M)<0 (X-p) +1 



eo(A) - h. (29) 



All of these steps are standard, with the exception of 
Eq. (24), leading to the equation Eq. (29) which is again 
standard with the exception of its nonstandard right- 
hand side. Here h is the chemical potential, and the 
integral is taken over only the region of //, where f (//,) < 0. 
Solving this equation for e(A) produces the excitation 
spectrum of the system, which is the quantity we would 
like to compute. Notice that the coupling constant g is 
not explicitly present, except through the definition of 
eo(A) in Eq. (14). 

If e(A) > 0, then the excitation at this A is a particle. 
If, on the other hand, e(A) < 0, then the excitation is a 
hole whose energy is — e(A). 



III. THE EXCITATION SPECTRUM 
A. Dimensionless parameters 

The excitation spectrum can be found by solving the 
equation Eq. (29). This can only be done numerically. 
To do this in a meaningful way, let us first study the 
scale of the parameters involved in Eq. (29). 

Since the interactions occur at a finite range ro , we will 
restrict the possible values of momenta k to the range 



k e 



(30) 



(as if the model Eq. (2) is defined on a lattice of lattice 
spacing tq). We would also like to make sure that the in- 
teractions are sufficiently weak so that particles moving 
with momenta close to n/ro would be close to noninter- 

acting. This can be achieved if the Ag"^ /ro is much smaller 
than {v — u)^k^ where k ~ 7r/ro in Eq. (10). This gives 



{u — vy 



<C 1. 



(31) 



From now on, we adopt this assumption. 

Second, it is convenient to rescale the rapidity A to 
simplify the expression for the energy spectrum eo(A). 
We introduce a parameter 



25v^ 
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< 1 



(32) 
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and define 



A = Ac. 



(33) 



We also introduce the dimensionless rescaled energy spec- 
trum 

_ \v-u\^ _ \v-u\^ 
e = ^ — e, eo = — — cq. (34) 



gy/uv gvuv 
The equation Eq. (29) gets simplified to 

eCX)-if d^-^^^^^^=X-l-h. (35) 
Je{il)<0 (A - /i)2 + c2 A 

This equation can only be solved numerically, even in the 
physical limit of c <C 1. 

The range of allowed momenta provides a natural cut- 
off for A. Eq. (12) together with Eq. (30) gives 
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In turn, this gives for A 

Ae 
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The integration range over fb in Eq. (35) is over these two 
combined intervals. 

Now wc arc in the position to solve Eq. (35) numeri- 
cally. The standard method is by interacting the relation 
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This leads to e„(A) quickly diverging to negative infin- 
ity as n increases. And indeed, the proof given in Ref. 
[11] regarding the convergence of this procedure is not 
applicable to Eq. (35). 

Instead, we use a different technique. We define a func- 
tional 
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is equivalent to Eq. (35), up to a multiplication by 
1 — signe(A). Then we introduce an extra fictitious pa- 
rameter T, and construct the solution to the equation 



di{X,T) 
dr 



SQ 
'6e{X,T) 
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in the limit where r — >■ oo. 

This procedure allows us to compute e(A) for all such 
A that e(A) < 0. One drawback of this procedure is 
that once e(A,T) = for some A and some r, it will 
remain zero for larger r. As a result, ?(A) can become 
"trapped" at zero whereas it might actually be negative. 
We fix this problem by supplementing it with iterations 
(38). Once the initial e(A) used for iterations is close to 
the solution of Eq. (35), subsequent iterations will not 
diverge. Indeed, suppose 



e(A)=e,(A) + 5e(A), 



(42) 



where is the solution of Eq. (35), and where Si <^ Cg. 
Then we find 



den+i{X) = 
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It is now fairly straightforward to prove that 

r dX (Sin+iWY < r dX (6~en{X)f , (44) 

which proves that the iteration procedure is not diver- 
gent. 

Once we do that, we construct the rest of this function 
by using Eq. (35) as a definition of e(A), or 

e(A) = i^/. d, (45) 
Je{X)<0 (A - /i)2 -I- c2 A 

B. Numerical solution 

We now use this procedure to construct solutions to 
Eq. (35). We take representative parameter values 



c = 0.1, - =2. 

u 



We then take initial value 



e(A,r) 



-1. 



r=0 



(46) 



(47) 



Then we run the Eq. (41) in steps of dr = 0.01 up to 
T = 100. The integrals are computed by discretizing the 
range of A into 3200 intervals. After that, we use e(A, 100) 
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as an input to the iteration procedure Eq. (38) where we 
iterate only once. This seems to be enough to generate a 
solution of Eq. (35) with a reasonable accuracy of about 
10^"^. The accuracy is defined as 




FIG. 1: ?(A) for c = 
( 


= 0.1, vju = 2, h = 0. 
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FIG. 2: Same figure as in Fig. 1 but the region close to A = 
enlarged. 

First we illustrate the solution for h — 0. Fig. 1 shows 
e(A) for the entire range of A. Fig. 2 shows the region of 
small A enlarged. 

One sees that at A < 0, e < 0. So for the branch of exci- 
tation spectrum at negative rapidities, all the excitations 
are holes, and they are all gapped. This is illustrated on 
Fig. 3. For A > 0, e changes sign at some value of A. So 
here we have gapless right moving excitations with the 
linear spectrum. This is illustrated on Fig. 4. 

By decreasing the chemical potential h, it is possible 
to make the A < molecules massless, whiling keeping 
the A > molecules massless as well. This is illustrated 
on Figs. 5 and 6. 

This last observation implies that the system we study 
undergoes a phase transition (or possibly more than one) 



e 




FIG. 3: Same figure as in Fig. 1 but the region of A < close 
to A = enlarged. 



e 




FIG. 4: Same figure as in Fig. 1 but the region of A > close 
to e = enlarged. 

as a function of the chemical potential (or of density), 
with the massless mode acquiring a gap. This phase 
transition looks similar to the transition observed in the 
problem defined by Eq. (1) in Ref. [7]. To elucidate the 
nature of this transition, it seems useful to study Eq. (2) 
using the bozonization techniques which we leave as a 
subject for future work. 
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FIG. 5: Same figure as in Fig. 3 but now the chemical poten- 
tial h = -50. 
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FIG. 6: Same figure as in Fig. 4 but now tiie chemical poten- 
tial h = -50. 



IV. CONCLUSIONS 

We have developed an exact solution to the problem 
of chiral atoms and molecules propagating in one dimen- 
sions with interactions controlled by a Feshbach reso- 
nance. The remaining outstanding issue is whether this 
Hamiltonian can be realized in a realistic cold atom ex- 



periment. 

A particular tantalizing question is whether the 
fermionic analog of Eq. (2) (the one where and d are 
fermionic creation and annihilation operators, and where 
an additional "spin" index is necessary to make the qubic 
term non-zero) is integrable. If so, the development of an 
exact solution of such a problem would be an interesting 
direction of further research. The fermionic version of the 
3- wave model was shown to be integrable in Ref. [12]. 
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